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ABSTRACT 



Context. According to the modern cosmological paradigm, galaxies and galaxy systems form from tiny density perturbations gener- 
ated during the very early phase of the evolution of the Universe. 

Aims. Using numerical simulations, we study the evolution of the density perturbation phases of different scales to understand the 
formation and evolution of the cosmic web. 

Methods. We apply the wavelet analysis to follow the evolution of high-density regions (clusters and superclusters) of the cosmic 
web. 

Results. We show that the maxima and minima positions of density waves (their spatial phases) almost do not change during the 
evolution of the structure. Positions of density perturbation extrema of are more stable for large scale perturbations. In the context 
of the present study we call density waves of scale > 64 h~ l Mpc large, waves of scale ^ 32 h~ l Mpc medium, and waves of scale 

8 h~ l Mpc small, within a factor of 2. 
Conclusions. In the cosmic structure formation of the synchronisation (coupling) of density waves of different scales plays an impor- 
tant role. 

Key words, cosmology: large-scale structure of the Universe; cosmology: early Universe; cosmology: theory; methods numerical 
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1. Introduction 

The basic structural elements of the Universe are filamentary 
superclusters and voids f orming a web-like structure - the 
supercluster-void network (lEinasto et al.lll980(: | Zeldovich et al.l 
ll982Hde Lapparent et al.lll986HBond et al.lll996l) . The standard 
I cosmological paradigm pr edicts that a period of accelera t ed ex- 
pansion, dubbed inflation dStaro binskv 1 98(1 lGuth|[l98 li iLindd 
ll982UAlbrecht & Steinhar dt 1982), generated density fluctua- 
tions (iMukhanov & Chibisov1ll98ll:[Hawkinglll982tlStarobinskvl 
ll982tlGuth & Pill 1982b as well as primordial gravitational waves 
dStarobin skv 1979) through quantum-gravitational processes. In 
the simplest form of this scenario, the primordial density field 
is predicted to form a statistically homogeneous, isotropic and 
almost-Gaussian random field, after the transition from quantum 
to approximate classical description of perturbations. 

If the hypothesis of primordial Gaussianity is correct, then 
density waves of different scales began with random and uncor- 
rected spatial phases. As the density waves evolve, they interact 
with others in a non-linear way. This interaction leads to the gen- 
eration of non-random and correlated phases which form a spa- 
tial pattern of the present cosmic web. Owing to non-linear pro- 
cesses during galaxy formation and the physical biasing prob- 
lem (almost no galaxies form in low-density regions), the present 
density field is highly non-Gaussian. There have been a number 
of attempts to gain quantitative information on the behaviour of 
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phases in gravitational systems; for a review see IColesI ( 2009) 
and references therein. 

Using numerical simulations, iRyden & Gramannl (Il99ll) 
showed that initial phase inf ormation is rapidly lost in short 
wavelengths during evolution. Hi kage et al.l (120051) analysed the 
clustering of SDSS galaxies using the distribution function of 
the sum of Fourier phases. Fourier phases are statistically inde- 
pendent of the Fourier amplitudes, thus the phase statistics plays 
a complementary role to the conventional two-point statistics of 
galaxy clustering. Gaussian fields have a uniform distribution of 
the Fourier phases over < 6k < 2n. Therefore, characteris- 
ing the correlation of phases is expected to be a direct means to 
explore non-Gaussian features. From the comparison of obser- 
vations with mock catalogues constructed from N-body simula- 
tions, the authors find that the observed phase correlations for 
the galaxies agree well with those predicted by the spatially flat 
ACDM model, evolved from Gaussian initial conditions. 

The analysis in Fourier space is, however, not sensitive to 
the location of particular high-density features in real space, 
such as filaments, clusters, and superclusters. To have a bet- 
ter understanding of the texture of the cosmic web, the web 
must be studied in the real space. Different statistical measures 
have been used to d escribe quantitatively the cosmic te xture, 
for recent reviews seelMartmez & Saarf (l2002h . lSaarf (120091) . and 
Ivan de Wevgaert & Schaapl (120091) . 

One of these statistics is the wavelet analysis, which a naly- 
ses p roperties of waves of various scales in real space (see lJonesI 
2009, and references therein). Wavelet analysis has been used to 
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detect voids and filaments in the C enter for Astrophysics (CfA) 
survey fir st slice ( Slezak et "aDl993h . to de-noise the galaxy dis- 
tribution dMartmez et al.ll2QQ5b . to detect the integrated Sachs- 
Wolfe e ffect in the cosm ic microwave background (CMB) ra- 
diation dVielya et ani2Q06l), to s tudy the discreteness effects in 
simulations dRomeo et al.ll2008l) . and for many other purposes 
where the spatial position of structural elements is important. 

The goal of this paper is to investigate the evolution of 
the texture of the cosmic web. It is well known that on small 
scales the original phase infor mation is lost durin g the non-linear 
stage of the evolution (iRyden & Gramannlll99ll) . On the other 
hand, the main large-scale skeleton of the texture of the cosmic 
web i s already determined by t he initial gravitational potential 
field (iKofman & Shandarinlll988l) . If considered from the point 
of view of the Fourier decomposition, maxima and minima in 
any spatial distribution occur in the points where phases of the 
Fourier modes are synchronised. However, there is a non-trivial 
problem, which is not completely solved yet. 

In the classical "down-top" model, like the isocurvature 
model, there are no built-in large-scale features. The structure 
formation starts from small-scale systems, which grow by ran- 
dom clustering. In the classical "top-down" model, like the adia- 
batic neutrino-dominated hot dark matter model, there is a built- 
in cut-off scale, which determines the scale of the structure. The 
presently accepted dark energy dominated ACDM model is es- 
sentially a "down-top" model, because the structure formation 
starts from small systems. This model differs from the classi- 
cal "down-top" model in one important detail - here a broad 
power spectrum of density perturbations is present. Objects of a 
smaller scale and mass form earlier. But in the case of a broad 
power spectrum of perturbations, it becomes a non-trivial ques- 
tion whether extrema of perturbations of a given scale will re- 
main at the same places if perturbations of larger scales are 
added. This may occur only if some phase synchronisation or 
coupling between perturbations of different scales exists. In the 
case of the broad wave spectrum, extrema of density perturba- 
tions should define locations where gravitationally bound ob- 
jects and voids form first. On the other hand, the gravitational 
potential defines the location of the skeleton of the cosmic web 
knots. Consequently, it is not clear at all why extrema of den- 
sity perturbations coincide with knots defined by the gravita- 
tional potential. In other words: Why is the skeleton stable in 
the "down-top" ACDM model? As we see in this paper, just be- 
cause of this synchronisation between waves of different scales. 

Studies of Fourier phases show that the phase coupling in 
the non-linear regime plays an importa nt role in the form ation 
of the fine texture of the cosmic web (iChiang et al .1120021) . To 
avoid complications caused by the highly non-linear regime, we 
concentrate on the evolution of waves at intermediate and large 
scales using the wavelet decomposition of the evolving density 
field. The Fourier modes are fully specified by their wavelength, 
their orientation, and phase. Because the phase determines where 
the maxima and minima are located along a Fourier mode, we 
also use the same terminology (somewhat less strict this time) 
once we speak about the wavelet decomposition of the density 
field. Here we assume that the (spatial) phase and the locations 
of the maxima and minima carry the same information, and thus 
will use these terms interchangeably in the following. Also, quite 
often the locations of the cells inside the cubical density grid are 
located as (/, j, k), with i, j, and k are integers that run from 1 to 
N, where N - 256 is most often assumed throughout the work. 

We shall focus our attention on the two main problems: the 
evolution of phases (positions of maxima) of density perturba- 
tions at medium and large scales, and the phase coupling (syn- 



chronisation) of perturbations of different scales. To follow the 
evolution of perturbations of different size, we use simulations 
in boxes of various sizes from 100 to 768 h~ l Mpc. To find 
the sensitivity of our results to the resolution, we make sim- 
ulations with 25 6 3 and 512 3 cells and equal number of parti- 
cles. For comparison with the real Universe, we shall calcu- 
late the density field and its wavelet decompositions for a slice 
(wedge) of the Sloan Digital Sky Survey (SDSS). Preliminary 
results of this s tu dy ha ve been reported at several conferences 
(Einasto 2006a : b, 2009). This paper is a follow-up of a study by 
Ein asto et al.l (l2005¥ of the environmental effects of clusters in 
SDSS and simulations. 

In the next section we describe the numerical models used 
in this study. We also make a qualitative wavelet analysis of 
the simulated density field, follow the density evolution in time, 
and compare the evolution of density waves of various scales. In 
Section 3 we make a correlation analysis of wavelet-decomposed 
density fields. In Section 4 we analyse the luminosity density 
field of the SDSS, and study the role of phase synchronisation 
in the formation of clusters and superclusters. In Section 5 we 
discuss our results. The last section gives our main conclusions. 



2. Qualitative wavelet analysis of the cosmic web 
evolution 

2.1. Modelling the cosmic web evolution 

In order to understand the evolution of the supercluster-void 
phenomenon, numerical simulations need to be performed in a 
box that contains both small- and large-scale waves. The most 
common systems of galaxies are groups and galaxy clusters 
with a characteristic scale of ~ 1 h~ l Mpc, therefore the sim- 
ulation must have at least a resolution of this scale. On the 
other hand, the largest non-percolating systems are superclus- 
ters, which have a characteristic scale up to ~ 100 h~ l Mpc. 
Superclusters have a very different richness from small systems 
like the Local supercluster to ver y rich systems like the Shapley 
supercluster (lEinasto et al .1120011) . Clearly this variety has its ori- 
gin in density perturbations of still larger scales. Thus, to under- 
stand the supercluster-void phenomenon correctly, the influence 
of very large-scale density perturbations should be studied, too. 

To have both high spatial resolution and the presence of 
density perturbations in a broad interval of scales, we used 
in this analysis simulations in boxes of sizes 100 h~ l Mpc, 
256 h~ l Mpc, and 768 h~ l Mpc, with the resolution of 256 3 
particles and simulation cells; these models are designated as 
Ml 00, M256 and M768, respectively . For simulations we used 
the AMIGA code dKnebe et al.l 120011) . This code uses an adap- 
tive mesh technique in regions where the density exceeds a 
fixed threshold. For comparison we used a model of box size 
256 h~ l Mpc with the resolution of 512 3 particles and cells, this 
model is designated L2 56. In the last ca se we used the GADGET 
code for simulations (Spr ingell 120051) . Results obtained from 
models with different resolution are very similar in the context of 
the present study, thus in the graphical representation of our re- 
sults we mostly use models Ml 00, M256, and M768. The model 
L256 has been used to find a catalogue of density field (DF) clus- 
ters of galaxies, to study the dependence of the mass of DF clus- 
ters on the density of the environment, and to compare wavelets 
of models with different cut-off scale (see Figure [8] below). 

The initial density fluct uation spectrum wa s generated us- 
ing the COSMICS code bv iBertschingerl (I1995IPI We assumed 
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Fig. 1. High-resolution density fields and their wavelet decompositions for models Ml 00 and M768, upper two and lower two rows, 
respectively. Model Ml 00 fields are given at k = 1 coordinate, M768 fields at k = 217 coordinate. Left panels show the density 
fields, next panels their wavelet decompositions. For both models upper panels show the field and wavelets at redshift z = 0, and 
lower panels at redshift z= 100. For model Ml 00 we show the wavelets w6, w5, and w4, for model M768 wavelets w5, w4, and w3. 
In plotting density fields and wavelets we use only over-dense regions: i.e., lower density limit (which corresponds to black level 
in DS9 colour palette SLS) is taken equal to 1 for the density field, and for wavelets. Upper density levels (which correspond to 
white in palette SLS) are given in Table[TJ Table[T]gives also the upper smoothing wavelength Sw^d in h~ l Mpc for all wavelets. 



cosmological parameters Q m = 0.28, Da = 0.72, erg = 0.84, 
and the dimensionless Hubble constant h = 0.71; to generate 
the initial data we took the baryonic matter density £2b = 0.044. 
Calculations started in an early epoch, z = 500. Particle posi- 
tions and velocities were extracted for 13 epochs between red- 
shifts z— 100, . . . , 0. In the present study we used only a part of 
these data, depending on the goal of the task. 



2.2. Wavelet decomposition of simulated density fields 

To see the effects of density waves of different scales and to 
understand the evolution of the density field, we shall use our 
models M100, M256, and M768. The density fields were cal- 
culated for all simulation epochs and were used to find the 
wavelets up to the level 6. The analysis was made in three di- 
mensions. For illustrations of the results we use two-dimensional 
slices in i, j-coordinates through the simulation box at fixed 
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Fig. 2. High-resolution density field of the model M256 is shown in the left column, atk= 153 coordinate. The second, third, and 
fourth columns shows the wavelet w6, w5, and w4 decompositions at the same k, respectively. The upper row gives data for present 
epoch, z = 0, the second row for redshift z = 1, the third row for redshift z = 5, and the last row for redshift z = 10. Densities are 
expressed in linear scale. Parameters used for plotting are given in Table [T] 



^-coordinates. The calculation of wavelets is explained in 
Appendix A. 3. 

The evolution of waves in models Ml 00 and M768 are 
shown in Fig.[T] High-resolution density fields at redshifts z = 
and z = 100 are shown in the left panels of Fig.[T] The wavelet 
decompositions of these density fields are shown for the same 
slices in the following panels of the figure. We show three 
wavelet levels: w6, w5, and w4 for the model Ml 00; and w5, 
w4, and w3 for the model M768. These wavelets show the evo- 
lution of large and intermediate scale waves. Upper smoothing 
wavelengths and upper density levels for plotting are given in 
Table [TJ We note that the upper smoothing scale of a wavelet 



wn is equal to the kernel diameter, used in the calculation of the 
density field D n , see Appendix A. 3. 

To see the evolution of the density field at intermediate time- 
steps we show in Fig.[2]the high-resolution density fields of the 
model M256 at four redshifts: z = 10, 5, 1, 0. Wavelet de- 
compositions at levels w6, w5, and w4 for the same redshifts are 
shown in the same Figure. Colour-coding of wavelets at different 
redshifts is chosen so that a certain colour approximately corre- 
sponds to the density level, corrected by the linear growth factor 
for that redshift. Parameters used for plotting are given in the 
Tableffl 
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Table 1. Parameters used in wavelet figures. 



Model 


z 


DFu 


Dw6u 


Dw5u 


Dw4u 


Swwu 


(1) 


(2) 


(3) 


(4) 


(5) 


(6) 


(7) 


M100 





40 


1.0 


1.8 


5.0 


50 


M100 


100 


1.15 


0.0085 


0.018 


0.03 


50 


M256 





20 


0.55 


1.5 


3.5 


128 


M256 


1 


10 


0.3 


0.75 


1.75 


128 


M256 


5 


4 


0.09 


0.2 


0.45 


128 


M256 


10 


2 


0.05 


0.1 


0.2 


128 


M768 





10 


0.25 


0.86 


2.0 


192 


M768 


100 


1.08 


0.004 


0.009 


0.018 


192 



Notes. 

1: Model; 
2: Redshift z\ 

3: Upper density limit of the high-resolution density field used in Figs[T] 

m 

4: Upper density limit of the wavelet w6 (w5 for the model M768) used 
in FigsIUEl 

5: Upper density limit of the wavelet w5 (w4) used in FigsHJO 

6: Upper density limit of the wavelet w4 (w3) used in FigsHJO 

7: Upper smoothing wavelength in h~ l Mpc of the largest wavelet Sw6u 

(Sw5u); the upper smoothing wavelength of the next wavelet Sw5u 

(Sw4u) is 2 times smaller than that of Sw6u (Sw5u), and that of Sw4u 

(Sw3u) is 4 times smaller than that of Sw6u (Sw5u). Lower smoothing 

wavelengths of all the wavelets are 2 times smaller than the upper ones. 



2.3. The evolution of density waves in time 

Let us first study the evolution of wavelets of various scales in 
time. It is well known that during the early stage of the evolution 
of the Universe the main evolution is in the growth of amplitudes 
of density perturbations. In the early stage of the evolution, the 
growth of amplitudes is almost linear. 

How well this approximation works in our numerical sim- 
ulations can be seen when we compare wavelets of different 
scales at various redshifts. The largest wavelets shown in Figs.Q] 
and [2] are w6 for the models M100 and M256, and w5 for the 
model M768. Upper smoothing scales of these largest wavelets 
are 50, 128, and 192 h~ l Mpc (models M100, M256, and M768, 
see Table [T]). The upper smoothing scales of next level wavelets 
for these models are 25, 64, 96 h~ l Mpc. The upper smoothing 
scales of the lowest level wavelets, shown in Figs. Q] and [2] for 
models M100, M256, and M768, are 12.5, 32, and 48 h~ l Mpc, 
respectively. The upper smoothing scales of the wavelet wl are 
1.56, 4, and 12 h~ l Mpc for these models. 

A look at Figs. [T]and [2] shows that in the model Ml 00 there 
are already considerable changes in positions and shapes of den- 
sity configurations in the wavelet w5 (upper smoothing scale 
25 h~ l Mpc). In the model M256 changes of the wavelet w5 (up- 
per smoothing scale 64 h~ l Mpc) are still relatively small. 

In the lowest scale wavelets (w4 of the model M256) the 
patterns of density distributions at various epochs are still fairly 
similar, however, changes in positions and density levels are 
more visible. These scales are upper limits of wavelengths of 
given wavelets, lower limits are twice smaller. Thus the mean 
smoothing scale of the wavelet w4 of the model M256 is 
32/ V2 = 22.6 /T 1 Mpc. Changes of patterns in the wavelet w4 
of the model Ml 00 (mean smoothing scale 8.8 /z _1 Mpc) are al- 
ready fairly strong, both in the position and the density level. 




Fig. 3. Evolution of the high-resolution density field of the model 
M256, at k = 153 coordinate. Upper panels are for z = and z = 
1, lower panels for z - 5 and z = 10. The densities are expressed 
in the logarithmic scale; this allows to see better the evolution in 
low-density regions. Only over-densities are shown, i.e., lower 
density limit for plotting is taken 1 . Upper densities for plotting 
with the DS9 package, corresponding to white, are 50, 25, 5, and 
2.5 for redshifts 0, 1, 5, and 10, respectively. Colour codes shown 
at the bottom correspond to redshift z = (upper left panel). 



To see the evolution of the filamentary density field better, 
we show in Fig. [3] the high-resolution density field of the model 
M256 again, with densities expressed in a logarithmic scale. We 
use the same epochs and k = 153 coordinate as in Fig. [2] In 
the lower left region a very rich supercluster is located, visible 
in the wavelet w6 of Fig. [2] as a large white area. The position 
of the supercluster does not change with time, the density con- 
trast increases approximately proportional to the linear growth 
factor (used in setting colours for plotting). Most visible effects 
are the contraction of the filamentary system towards the centre 
of the supercluster, thinning of filaments, and merging of small 
filaments. 

Our simple qualitative analysis shows that the evolution of 
waves of medium and large scales is approximately linear as 
expected. Positions of maxima of these waves change very lit- 
tle with time. Waves of smaller scales change their positions 
of maxima, this can be considered as phase shifts of maxima. 
The shift is larger for small-scale waves. These shifts are caused 
by two effects: the contraction of high-density regions (super- 
clusters), and the redistribution of matter on small scales, as 
illustrated in Fig. [3] The c ontraction of the h igh-density re- 
gions is visible in Fig. 17 by Tem pel etaD (120091) . which shows 
the changes of the sizes of particle samples, giving rise to the 
present-day clusters of galaxies. In central regions of superclus- 
ters clusters formed by merging of a large number of primordial 
clusters, collected from an approximately spherical volume of 
outer radius about 8 h~ l Mpc. The lower the global density, the 
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wavelet level 

Fig. 4. Behaviour of the correlation coefficient for different red- 
shift pairs (Zji, zj) e {(30, 10); (30, 5); (30, 1); (30, 0)} for all six 
wavelet levels. We see that for the largest smoothing scale, i.e. 
w6, all the correlators stay very close to r = 1 . At lower wavelet 
levels, as the redshift drops below z = 30, the lines start to de- 
viate from r = 1. Thus, on the largest scales the information is 
approximately preserved, while on the smallest scales the infor- 
mation is gradually erased. 



z=5 — * — 
z=10 x 




1 2 3 4 5 6 

wavelet level 



Fig. 5. Correlators at fixed redshifts z\ - Zj = 0, 1, 5, 10, 30 from 
top to bottom (from red to yellow) for high-density regions. Note 
that in this figure we define the high-density regions as the ones 
where the pixel values on the largest smoothed wavelet field 
reach the top 10% of all the values. Using these pixels, we de- 
fine the mask, which is used for correlating only the highest den- 
sity pixels with the corresponding values smoothed on a smaller 
scale. 



smaller is the radius of the sphere from which present clusters 
are formed. 



2.4. The comparison of evolution of density waves of various 
wavelengths 

Now we compare the evolution of density waves of different 
wavelengths, using the wavelet decomposition of our model den- 
sity fields at various epochs. 

A look at Fig. [2] shows that at all redshifts high-density 
peaks of wavelets of medium and large scales almost coincide. 
In wavelets of smaller scales, there is sometimes more than one 
peak in the high-density region of the next larger scale, but at 
least one peak is present there in practically all cases. In other 
words, wavelets of various scales have a tendency of phase cou- 
pling or synchronisation in peak positions. Positions of high- 
density peaks of wavelets wl and w2 coincide for any fixed red- 
shifts (not shown in Figs. Q] and [2l but seen in original wavelet 
figures). The highest of these small size density peaks also co- 
incides in position with peaks of higher level wavelets. We con- 
clude that the growth of small-scale peaks is amplified in high- 
density regions of large waves: peaks in w4 to w6 in Fig. [2] and 
peaks of wavelets of various scale in Fig. 7 of lEinasto et al.l 
( 201 H). This amplification of density peaks of perturbations of 
various scales near peak positions is the reason for wave cou- 
pling (synchronisation). 

The general conclusion from the wavelet analysis of all 
scales and from the comparison with corresponding density 
fields is that the synchronisation of peak positions of wavelets 
of various scales represents a general property of the evolution 
of the density field of the Universe. A quantitative analysis of 
this effect is given in the next section. 



3. Correlation analysis of wavelet-decomposed 
density fields 

In this section we attempt to quantify some of the qualitative 
statements given above. To this end we perform the correlation 
analysis of wavelet-decomposed density fields. For clarity we 
will use only the model M256 because the other ones lead to very 
similar results. Here we have six wavelet levels: wl, w2, . . . , w6 
with the upper smoothing scales of 4, 8, . . ., 128 h~ l Mpc, re- 
spectively. We use the simulated density fields at five different 
redshifts: z = 30, 10,5, 1,0. 

Most generally, we have chosen two redshifts, Zi and Zj, 
along with two wavelet levels, w m and w n , and calculated the 
following correlators: 

r w m Zi,w n Zj — i • \I) 

Here 6 WmZ . corresponds to the wavelet-decomposed density field 
for level w m at redshift Zi. The angle brackets represent ensem- 
ble average, which under the ergodicity assumption is replaced 
by a simple spatial average. Note that all wavelet-filtered fields 
have zero mean, i.e. (S WmZi ) = 0. Also note that we calculate the 
correlators for zero lag only, i.e., we do not shift one field with 
respect to the other one. In the following we will use two types 
of correlators: 

1. Fixed wavelet scale correlators, i.e., w m - w n , at different 
redshifts. 

2. Correlators at fixed redshifts, i.e., Zi = Zj, for different 
wavelet levels. 

In the first case we take one of the density fields always at 
redshift z - 30, which is high enough for all the scales of interest 
to be well in the linear regime. It is easy to understand that under 
the linear evolution, where the values of the density contrast 6 get 
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just multiplied by the same scale-independent but time- varying 
factor, the correlation coefficient should always stay at the value 
r = 1, i.e., all the initial information is well preserved. In Fig. [4] 
we show the behaviour of the correlation coefficient for different 
redshift pairs: (zji, zj) e {(30, 10); (30, 5); (30, 1); (30, 0)} for all 
the six wavelet levels. It is easy to see that for the largest smooth- 
ing scale, i.e., w6, all the correlators stay very close to r = 1, 
while later on, as the other redshift drops below than z = 30, 
the lines start to decline from r = 1, especially at the smallest 
scales. Thus, on the largest scales the information is approxi- 
mately conserved, while on the smallest scales the information 
gets erased. The lower the redshift of the other density field, the 
higher the effect of the information erasure. In practice, for the 
cases z = 10 and z = 5 the loss of information is relatively mod- 
est for all wavelet levels. For z = 1 and z = the information is 
approximately preserved if the wavelet level w > 4. 

In Fig. \5\ we plot the correlators at fixed redshifts Zi = 
Zj = 0, 1, 5, 10, 30 from top to bottom (from red to yellow). For 
w m = 1 the curves are peaked at w n = 1, while they gradually 
drop as w n is increased to higher values. Similarly, the curves 
for w m = 2 are peaked at w n = 2, and get reduced as the dis- 
tance increases from this point. For the other w n values the be- 
haviour is very similar. As long as the evolution proceeds in a 
linear manner, i.e., the growth depends only on redshift, but is 
independent of the wavelet scale, the coupling kernels plotted in 
Fig. \5\ should stay exactly the same. However, we see that the 
lower the output redshift, the broader are the coupling kernels, 
i.e., a nonlinear evolution in general leads to the additional cou- 
pling of the nearby wavelet modes. Also, note that in this figure 
we define the high-density regions as the ones where the pixel 
values on the largest scale wavelet field reach to the top 10% 
of the values. Using these pixels we define the mask, which is 
used for correlating only the highest 10% of the density values 
with the corresponding values smoothed density on the smaller 
scale. Thus, in this case we have looked only how the highest 
density regions correlate with each other. We note that owing to 
nonlinear evolution, the coupling kernels become broader and 
broader with time. It is important to realise that even for the lin- 
ear evolution of the Gaussian density field the nearby wavelet 
levels at fixed redshift are significantly coupled, because in this 
case the neighbouring levels tend to contain some of the com- 
mon Fourier space modes. Assuming only the linear evolution, 
clearly the coupling does not change with redshift. 

4. Wavelet analysis of the SDSS luminosity density 
field 

4.1. The SDSS luminosity density field 

The luminosity density field was calculated using galaxy data of 
the Main sample of th e contiguous Northern area of the SDSS 
data release 7 (DR7) (lAbazaiian et al.ll20Q9b . The DR7 sample 
has, after applying extinction corrections, Petrosian r-magnitude 
limits, 14.5 < r < 17.77. We used for this analysis galaxies in the 
redshift interval 2700 < cz < 60000 km s _1 . The SDSS data re- 
duction procedure consists of two steps: (1) the calculation of the 
distance, the absolute magnitude, and the weight factor for each 
galaxy of the sample, and (2) the calculation of the luminosity 
density field using an appropriate kernel and a chosen smoothing 
length. We estimated total luminosities of groups and isolated 
galaxies in our flux-limited sample, using weights of galaxies 
that take into account galaxies and galaxy groups too faint to 
fall into the observational window of absolute magnitudes at the 
distance of the galaxy. When calculating luminosities of galaxies 



we regard every galaxy as a visible member of a group within the 
visible range of absolute magnitudes. For details of the calcula- 
tion of the luminosity density field see Appendix A. 1 . A super- 
cluster catalogue ba sed on the luminosity d ensity field of SDSS 
DR7 is published by Liivama gi et al.1 (120101) . and the luminosity 
function of galaxies by Te mpel et al.l (1201 ll) . 
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Fig. 6. Upper left panel shows a two-dimensional rectangular re- 
gion of size 512x512 Mpc (using Hubble constant h = 0.8), ex- 
tracted from the high-resolution luminosity density field of the 
Northern equatorial wedge of the SDSS. The density was cal- 
culated using Gaussian smoothing with the rms scale 0.8 Mpc. 
The observer is located at the lower left corner. Colour-coded 
density levels used in plotting are in the interval from to 10 in 
mean density units with the highest value corresponding to the 
white colour, using the SAOImage program DS9 colour palette 
SLS. The upper right panel shows the same density field as in 
the previous panel, but the phases randomly shifted. With shift- 
ing the phases of density waves some densities become negative, 
thus in this case colour codes are in the density interval ±3.5. 
Middle and lower panels show the wavelets w7, w6, w5, and 
w4 of the wedge. The middle left panel shows density waves of 
scales approximately equal to half of the box size, each follow- 
ing panel shows waves with half the length. Wavelet densities 
used in plotting have both negative and positive values, i.e., they 
correspond to under- and over-densities of matter, colour cod- 
in g is linear. Superclu ster numbers according to the catalogue 
by Ein asto etaD (120031) are shown in the middle right panel w6. 
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For the present wavelet analysis we used the Northern equa- 
torial wedge of 2.5 degrees thickness up to redshift z = 0.2, see 
Fig. 00. The high-resolution projected two-dimensional luminos- 
ity density field was found using Gaussian smoothing with scale 
0.8 Mpc. Densities were found in a grid of bin size 1 Mpc. To 
reduce the wedge to unit thickness, density values were divided 
by the ratio of the thickness of the wedge at a particular distance 
from the observer, and at the mean distance. We use a rectan- 
gular sub-region of the field suitable for the Fourier and wavelet 
analysis. 

4.2. Fourier and wavelet decomposition of the SDSS density 
field 

The importance of the phase information in the cosmic web 
formation has b een understood long ago, as demonstrated by 
IColes & Chiangl (l2000h by randomising phases of a simulated 
filamentary network. To study the role of phase information in 
more detail, we extracted a rectangular region with the box size 
512 Mpc, calculated for the Hubble constant h = 0.8, from 
the density field of the SDSS Northern equatorial slice. This 
region is shown in the upper left panel of Fig. We Fourier- 
transformed the 2D density field and randomised phases of all 
Fourier components, and thereafter Fourier- transformed it back 
to see the resul ting density field . A sim ilar procedure has been 
applied also by lColes & Chiangl d2QQQh . The modified field has 
the same amplitudes of all wave-numbers k as the original field, 
only the phases of waves are different. Results are shown in the 
right upper panel of Fig. We see that the whole structure of su- 
perclusters, filaments and voids has gone, the field is fully cov- 
ered by tiny randomly spaced density enhancements. There are 
no clusters of galaxies in this picture, comparable in the luminos- 
ity to real clusters of galaxies. Fourier ph ase rand omising in nu- 
merica l models by lChiang & Colesl (l200Qh and lColes & Ch iang 
(2000) shows similar results. 

This simple example shows the importance of the density 
perturbation phases in the cosmic web formation. 

4.3. Phase synchronisation and the supercluster structure 

Next we use the wavelet analysis to investigate the role of den- 
sity waves of different scales. Figure [6] shows wavelets 7 to 4 
of the Northern rectangular region. These wavelets characterise 
waves of length about 256, 128, 64, and 32 Mpc, respectively 
(note that all scales in the Figure [6] are given using Hubble con- 
stant h = 0.8). In wavelet figures both under- and over-densities 
are shown. Extreme levels were chosen so that main features of 
the structure are well visible. 

The middle left panel of Fig. [6] shows the waves of length 
about 256 Mpc. In its highest density re gions there are three very 
rich superclusters: N20 from the list by lEinasto et al.l ([2003), lo- 
cated in the upp er part of the figure, supercluster N13 (SCL126 
from the list by lEinasto et aD (1200 ll) in the Sloan Great Wall) 
near the centre, and supercluster N02 (SCL82) in the lower right 
part of the panel. 

The next panel shows waves of scales about 128 Mpc. Here 
the most prominent features are superclusters N13 (SCL126) and 
N02 (SCL82), also the supercluster N23 (SCL155) in the upper 
left part of the panel is fairly strong, seen as a weak density peak 
already in the previous panel. In addition we see the supercluster 
N15 just above N13 near the minimum of the wave of 256 Mpc 



2 This wavelet decomposition was first presented by Einasto ( 2006b). 



scale, and a number of poorer superclusters located mostly in 
voids defined by waves of larger size. 

The lower left panel plots waves of scales about 64 Mpc. 
Here all superclusters seen on larger scales are also visible. A 
large fraction of density enhancements are either situated just in 
the middle of low-density regions of the previous panel, or they 
divide massive superclusters into smaller subunits. This property 
is repeated in the next panel. Here the highest peaks are substruc- 
tures of rich superclusters, and there are numerous smaller den- 
sity enhancements (clusters) between the peaks of the previous 
panel. 

When we compare the density waves of all scales, we ar- 
rive at the conclusion that superclusters form in regions where 
density waves of medium and large scales combine in similar 
over-density phases. The larger the scale of the wave where this 
coincidence takes place, the richer the supercluster. 

Similarly, voids form in regions where density waves of 
medium and large scales combine in similar under-density 
phases. In large voids medium-scale perturbations generate a 
web of filamentary structur es with knots; for a description of the 
formation of such a web see Bond et al. (1996). The influence of 
density perturbations of various scale s on the cosmic web fo r- 
mation is discussed in more detail by Suhh onenko et al.l (1201 ll) . 

This simple analysis very clearly demonstrates the role of 
phase coupling (synchronisation) of density waves of differ- 
ent scales in the formation of the supercluster- void network. A 
wavelet analysis of the full SDSS DR7 contiguous Northern re- 
gion is in progress, results will be published in a separate paper. 




Dglob 



Fig. 7. DF cluster masses of the model L256 as a function of the 
global density. Masses are given in solar units, densities in units 
of the mean density of the model. 



5. Discussion 

5.1. The role of phase synchronisation in the formation of 
clusters and superclusters 

The wavelet decomposition of SDSS data has shown that super- 
clusters form in regions where density waves of medium and 
large scales combine in similar over-density phases. Now we 
shall consider the role of phase synchronisation in the formation 
of galaxy clusters. 

Our previous analysis has shown that clusters of galaxies 
(haloes in simulations) form in places where small density waves 
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Fig. 8. Wavelets wl of models L256.256 and L256.008 at redshift z = 30 are shown in the left and right panel, respectively, at 
coordinate k = 1. Densities are expressed in linear scale, only over-density regions are shown. 



in a certain range of scales combine in similar over-density re- 
gions of waves. This is clearly seen as high level correlations of 
wavelets of levels 1, 2, and 3 in Fig. and close positions of 
maxima of wavelets w4 and w5, and positions of rich clusters in 
the high-resolution density field in Fig. [6) 

To see the role of phase synchronisation in the formation of 
clusters we shall use the density field (DF) clusters of the simu- 
lation of box length L = 256 h~ l Mpc with resolution 512 3 par- 
ticles and cells, the model L256. We define DF clusters as peaks 
of the high-resolution density field, calculated with the kernel of 
radius equal to the cell size of the simulation, the field Dq (see 
Appendix A. 3). DF cluster positions were identified as grid co- 
ordinates i, j, k of the local maxima of the field Dq. The mass of 
the DF clusters, M c \, was found by counting local density values 
in cells within ±4 cells from the central one, i.e., within a box 
of half-length 2.5 h~ l Mpc, centred on the peak. We express DF 
cluster masses in solar units using the masses of particles. T his 
sample of DF clusters was used by Suhho nenko et al.l (120111) to 
find DF cluster mass distribution and cluster-defined void radii. 
To investig ate the growth of structure s in the standard ACDM 
cosmogony jLudlow & Porcianil (|2010) used a slightly different 
algorithm to find density peaks in simulations. 

The relationship between DF cluster masses and global den- 
sities is shown in Fig. [71 Global densities at the peak positions 
were found using the density field calculated with a kernel of 
radius 8 h~ l Mpc, the field D4. The figure shows that there is 
an upper limit of masses of DF clusters for any given global 
density value. Global densities calculated wi th a kernel of radius 
8 /T 1 Mpc were used to define super clusters (lEinasto et al. 2003, 
I2006L 120071: iLiivamagi etafl l2Q10h . using a density threshold 
^ 4.5 in mean density units. We see that most massive DF clus- 
ters are located in superclusters. As the global density decreases, 
the upper mass limit of DF clusters also decreas es. A similar 
tende ncy is found in SDSS clusters of galaxies by lEinasto et al] 
( 2005): in high-density environment clusters have a mean lumi- 
nosity a factor of ^ 5 higher than in a low-density environment. 



This is the result of phase synchronisation: the larger the scale 
of density waves, where the maxima of waves of different scales 
have close positions, the larger are the masses of DF clusters. 
This also explains why there are no very rich clusters in a low- 
density environment: amplification by large-scale density waves 
is needed to form a very massive cluster. This amplification oc- 
curs only in the regions we call superclusters, as seen also in the 
SDSS data in Fig. 

Figure[7] shows that at any given global density value there is 
a large range of masses of DF clusters. This has a simple expla- 
nation. Small-mass DF clusters are located far from the maxima 
of large-scale waves, seen in wavelets of higher order. 

We conclude that in the formation of both small and large 
systems of galaxies the synchronisation of phases of density 
waves of different scales plays an important role. Our experi- 
ence shows that the wavelet analysis has a good diagnostic value 
in understanding the evolution of systems of galaxies of different 
richness. 



5.2. Phase synchronisation as a physicai process 

Phase synchronisation between different scales is a physical pro- 
cess that requires a causal connection between different points in 
space. Because standard inflation provides us with this connec- 
tion for scales up to very large ones, much exceeding the present 
Hubble scale, this synchronisation can already imprinted into the 
initial spectrum - actually, this is the standard modern paradigm. 
However, it has not yet been proved if this initial synchronisa- 
tion under the assumption of Gaussian statistics of perturbations 
in the linear regime is sufficient for the complete explanation of 
the structure and evolution of the cosmic web. Thus, our investi- 
gation may also shed light on the actual statistics of the primor- 
dial spectrum including its behaviour at large deviations from 
rms values. An alternative possibility might be some physical 
hydrodynamical processes acting before and during recombina- 



9 



J. Einasto et al.: Formation of the cosmic web 



tion. However, scales larger than the sound horizon at recombi- 
nation, ~ 146 Mpc accord ing to the most recent cosmological 
data by lJarosik et al.l (I2Q1QI) remain unaffected by any such pro- 
cess. 

Finally, larger scales may affect the evolution of smaller ones 
in the recent epoch. It is well known that density perturbations of 
large scales evolve almost linearly and do not change the i r posi- 
tions. This behaviour was used by lKofman & Shandarinl (Il988l) 
in the adhesion theory of the evolution of the large-scale struc- 
ture of the Universe. Applying Burgers equations, these authors 
predicted the actual geometrical structure of the cosmic web 
skeleton, using the primordial (post-inflation) gravitation poten- 
tial field. Numerical simulations with the same initial condi- 
tions confirmed the correctness of the prediction. This calcula- 
tion shows that the cosmic web skeleton was created at a very 
early stage of the evolution of the Universe. Our analysis has 
shown how this stability of the shape of the cosmic web skele- 
ton is expressed in wavelet terms. 

To understand better the synchronisation of waves of differ- 
ent scales as a physical process, we used two models of the 
L256 series. The first model has a full power spectrum, the 
model L256.256. The other model has a power spectrum trun- 
cated at wave-number k u so that the amplitude of the power 
spectrum on large scales is zero: P(k) = 0, if k < k u wave- 
length A t = 2n/k t = 8 /T 1 Mpc, the model L256.0Q8. Dat a of 
these simulations are described by Suhho nenko et al.l (1201 lh . 

The power spectrum of density perturbations has the highest 
power at small scales. Accordingly, the influence of small-scale 
perturbations relative to large-scale perturbations is strongest in 
the early period of structure evolution. For this reason the density 
fields and wavelets wl at early epochs are almost identical for 
the full model L256.256, and for the model cut at small scales, 
/i t = 8 h~ l Mpc, L256.008, see Fig. [8] Eventually, perturbations 
of larger scale start to affect the evolution. These perturbations 
amplify small-scale perturbations near maxima, and suppress 
small-scale perturbations near minima. In this way the growth 
of small-scale perturbations becomes non-linear. Thereafter still 
larger perturbations amplify smaller perturbation near their max- 
ima, and suppress smaller perturbations near their minima, and 
so on. The largest amplification (non-linearity) occurs in regions 
where maxima of perturbations of all scales happen to coincide. 
In such a way the synchronisation of phases of waves of different 
scales occurs as a natural process. The synchronisation of waves 
of different scales as a function of redshift z is seen graphically 
in Fig. 7 of lEinasto etaD (120111) . 



6. Conclusions 

Main conclusions of the present paper can be formulated as fol- 
lows. 

1. The wavelet analysis has demonstrated a good diagnostic 
value for studying the evolution of galaxy systems of vari- 
ous scales and masses. 

2. In the formation of cosmic structures the synchronisation 
(coupling) of density waves of different scales plays an im- 
portant role. 

3. Positions of density maxima of waves of large and medium 
scales practically do not change during the evolution. On 
smaller scales positions of density maxima change during 
the evolution, the changes are larger for waves with shorter 
wavelengths. 



4. Superclusters are objects where density waves of medium 
and large scales combine in similar phases to generate high- 
density regions. 

5. Voids are regions in space where density waves of medium 
and large scales combine in similar under-density phases. 

6. Clusters of galaxies are objects where density waves of small 
scales combine in similar over-density phases. 

7. The larger is the scale of the highest phase synchronisation, 
the richer are the clusters and superclusters. 
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Appendix A: Density field and wavelets 

A.1. Density field 

To estimate the expected total luminosity of groups or single 
galaxies, we assume that the luminosity functions derived for 
a representative volume can be applied also to individual groups 
and galaxies. Under this assumption the estimated total luminos- 
ity per one visible galaxy is 

hot = L ohs W d , (A.l) 

where L obs = L o 10 a4x(M O~ M) is the luminosity of the visible 
galaxy of absolute magnitude M (in units of the luminosity of 
the Sun, L©), and 

W d = ^- (A.2) 

£lf(Wl 

is the luminous-density weight (the ratio of the expected total 
luminosity to the expected luminosity in the visibility window). 
L\ and L2 are lower and upper limit of the luminosity window, 
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respectively. In our calculations we adopted the absolute mag- 
nitude of the Sun in the r filter M© = 4.64 (Blan ton & Roweisl 
l20Q7h . 

The ^-correction for the SDSS galaxies is calculated us- 
ing the KCORRECT alg orithm (version v4. 1 4) de veloped by 
IBlanton et aD (l2QQ3bh and lBlanton & Roweisl (120071). Evolution 
correcti on e has been calculated according to IBlanton et al.l 
(200 3^). For de tails of the data reduction procedure see 
iTago etal.1 (l2Q10h . 

In calculating of the total expected luminosity we used a 
double-power-law luminosity function with a smooth transition: 



0(L)dL oc (L/LT(1 + (L/L*yf- a)/7 d(L/L*), 



where a is the exponent at low luminosities (L/L*) <«c 1, 6 is the 
exponent at high luminosities (L/L*) » 1, y is a parameter that 
determines the speed of transition between the two power laws, 
and L* is the characteristic luminosity of the transition, simi- 
lar to the charact eristic luminosity of the Schechter function. As 
demonstrated by lTempel et al.l (120091) . the double-power-law lu- 
minosity function fits observed luminosity distribution of galax- 
ies and groups of galaxies better than the Schechter function, in 
particular in the high-luminosity end of the distribution. 

A.2. Kernel method 

To calculate a density field, we need to convert the spatial posi- 
tions of galaxies and their luminosities L t into spatial (lumi- 
nosity) densities. The standard approach is to use kernel densi- 
ties 



p(r) = V K(y-Yud)L u 



(A.4) 



where the sum is over all galaxies, and K (r; a) is a kernel func- 
tion of a width a. Good kernels for calculating densities on a 
spatial grid are generated by box splines Bj. Box splines are lo- 
cal and they are interpolating on a grid: 



(A.5) 



for any x and a small number of indices that give non-zero values 
for Bj(x). 

We use the popular B3 spline function: 



The three-dimensional kernel is given by the direct 
product of three one-dimensional kernels: 



Kf(r\ a, S) = K%\x\ a, S)K^\y; a, 6)K^\z\ a, 6), 



(A.8) 



where r ee {x,y,z}. Alth ough this is a direct product, it is 
isotropic to a good degree (Saar 2009). 

To calculate the high-resolution density field, we use the 
kernel of a scale equal to the cell size of the simulation, l c = 
Lb/N gr id, where Lb is the size of the simulation box, and Af gr id is 
the number of grid elements in one coordinate. 



(A.3) A.3. Wavelets 



We use the 'a trous' wavelet transform (Martinez & Saar 2002). 
The algorithm of the wavelet transform works as follows. Let 
us have a data set D (particles in simulations or luminosity 
weighted galaxies in SDSS data), located in a box of size nxnxn. 
The wavelet transform decomposes the data set as a superposi- 
tion of the form 



D 



^E 1 



(A.9) 



where Dj is the smoothed version of the original data D, and W j 
represents the structure o f D at scale 2 7 , (see lStarck et al.| [l998: 
IStarck & Mur tagh 2002). The wavelet decomposition output is 
/ three-dimensional density fields Dj and wavelets wj of size 
nxnxn. Following the traditional indexing convention, we de- 
note density fields and wavelets of the finest scale as j = 1. 
All density fields were calculated with the B3 spline kernel. The 
smoothed version of the original data, Dj = Dq, is the density 
field found with the kernel of the scale, equal to the cell size of 
the simulation l c . 

We calculated wavelets of index j by subtraction of density 
fields: 



D, 



(A. 10) 



where every higher level density field Dj was calculated with 
kernel size twice larger than the previous level field Dj-\. In such 
construction a wavelet of index j corresponds to density waves 
between scales A 7 _i = l c x 2 j ~ l/2 and Aj = l c x 2 j+l/2 , i.e., 
diameters of kernels used in calculation of density fields Dj-\ 
and Dj. 



Bs(x) = \z [\x - 2| 3 - 4\x - 1| 3 + 6|x| 3 - 4\x + 1| 3 + \x + 2| 3 ] .(A.6) Refer ences 
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We define the (one-dimensional) B3 box spline kernel of the 
width a as 



K^\x\a,d) = B 3 (x/a)(6/a), 



(A.7) 



where S is the grid step. This kernel differs from zero only in the 
interval x e [-2a, 2a\ \ it is close to a Gaussian with cr = 1 in 
the r egion x e [-a, a], so its effective width is 2a (see, e.g., Saar 
120091) . 

The kernel preserves the interpolation property exactly for 
all values of a and 5, where the ratio a/S is an integer. (This 
kernel can be used also if this ratio is not an integer, and a » 5\ 
the kernel sums to 1 in this case, too, with a very small error). 
This means that if we apply this kernel to N points on a one- 
dimensional grid, the sum of the densities over the grid is exactly 
N. 
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